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Abstract 



We use isomorphism ip between matrix algebras and simple orthogonal Clifford alge- 
bras C£(Q) to compute matrix exponential e A of a real, complex, and quaternionic 
^ — , , matrix A. The isomorphic image p — ip(A) in C£(Q), where the quadratic form Q 

has a suitable signature (p, q), is exponentiated modulo a minimal polynomial of p 
using Clifford exponential. Elements of C£(Q) are treated as symbolic multivariate 
polynomials in Grassmann monomials. Computations in C£(Q) are performed with a 
Maple package 'CLIFFORD'. Three examples of matrix exponentiation are given. 

-I— > 



1 Introduction 



Exponentiation of a numeric nxn matrix A is needed when solving a system of differential 
equations x' = Ax, x(0) = Xn, in order to represent its solution in a form e^*xo. It is well 
known that the exponential form of the solution remains valid when A is not diagonalizable, 
$— i ' provided the following definition of e A is adopted: 



eA = Z^TT' where A ° = L C 1 ) 

k=0 



n 



Equation ([j]) means that the sequence of partial sums S n = ^A fc /A;! — ► e A entrywise. 



fc=0 



Equivalently, (|l[) implies that \\S n — e A \\\ — > where \\A\\i denotes matrix 1-norm defined 
as the maximum of {||A/||i, j = l,...,n}, Aj is the jth column of a A, and ||-Aj||i is 
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n 

the 1-vector norm on C n defined as ||x||i = ^] \xj\. However, for several reasons, there 

i=l 

is no obvious wayQ to implement definition (||) on a computer, unless of course A is 
diagonalizable, that is, when A has a complete set of linearly independent eigenvectors 
(cf. 1). 

Another approach to solving x' = Ax is to find Jordan canonical form J of the matrix 
A. Let P be a nonsingular matrix such that P~ 1 AP = J. Then, if a change of basis is 
made such that x = Py, the matrix equation x' = Ax is transformed into y' = Jy and, at 
least theoretically, its solution is represented as e Jt c for some constant vector c. However, 
since the Jordan form is extremely discontinuous on a set of all n x n matrices, numeric 
computations of J are seriously ill-posed (cf. @, [|). 

In this paper we present another approach to exponentiate a matrix, let it be numeric 
or symbolic, with real, complex, or quaternionic entries, totally different from the linear 
algebra methods. It relies on the well-known isomorphism between matrix algebras over 
R, C, or H, and simple orthogonal Clifford algebras (cf. 0, ||, ||, 0]). This is not a matrix 
method in the sense that elements of the real Clifford algebra C£(Q) are not viewed here as 
matrices but instead they are treated as symbolic multivariate polynomials in some basis 
Grassmann monomials. This is possible due to the linear isomorphism Cl(V,Q) ~ f\V. 
The critical exponentiation is done in the real Clifford algebra C£ Ptq over Q with a suitable 
signature (p, q) depending whether the given matrix A has real, complex, or quaternionic 
entries. Three examples of computation of the matrix exponential with a Maple package 
'CLIFFORD' (cf. 1, 1, 0) are presented below. The Reader is encouraged to repeat 
these computations. 

In order to find matrix exponential e A , the following steps will be taken: 

- We will view elements of C£ Ptq as real multivariate polynomials in basis Grassmann 
or Clifford monomials. 

- We will find explicit spinor (left-regular) representation 7 of C£ p , q in a minimal left 
ideal S = C£ PA f generated by a primitive idempotent /. 

- For a matrix A (numeric or symbolic) in the matrix ring R(n), C(n) or H(n) where 
n = 2 m_1 , m = [^(p + q)] , we will find its isomorphic image p = p(A) in C£ Ptq !^ 

- We will find a real minimal polynomial p(x) of p and then a formal power series 
exp(p) modp(x) in Ci p ^ q . 

- We will check the truncation error of the power series exp(p) in C£ p ^ q via a polynomial 
norm, or in a matrix norm, both built into Maple^J 

- We will map exp(j?) back to the matrix ring M(n), C(n) or M(n) to get exp(^4). 

Before we proceed, let's recall certain useful facts about orthogonal Clifford algebras 
C£ P! q. For more information see Q. 

x It is possible to compute the exponential e At with a help of the Laplace transform method applied to 
an appropriate system of differential equations 
2 The brackets [ • ] denote the floor function 

3 It is also possible to use the R n topology where n' — 2 n , n — p + q. 
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- If p — q ^ 1 mod 4 then C£ Ptq is a simple algebra of dimension 2 n , n = p + q, 
isomorphic with a full matrix algebra with entries in R, C, or HI. 

- If p — q = 1 mod 4 then C^p,<j is a semi-simple algebra of dimension 2 n , n = p + q, 
containing two copies of a full matrix algebra with entries in R or H projected out 
by two central idempotents |(1 ± eie2 • • • e n ).[] 

- Ct P) q has a faithful representation as a matrix algebra with entries in R, C, H or 
R © R, HI © HI depending whether C£ Pj(? is simple or semisimple. 

- Any primitive idempotent / in C£ p>q is expressible as a product 

/ = i(l ± e Tl )\{l ± e T2 ) • • • i(l ± e Tfe ) (2) 

where {eji, er 2 , . . . , e^}, k = q — r q ^ p , is a set of commuting basis monomials with 
square 1, and r{ is the Radon-Hurwitz number defined by the recursion rj + 8 = r i + 4 
and 



i 





l 


2 


CO 
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5 


6 


7 


n 
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2 


2 


3 


3 


3 


3 



_ Ci P) q has a complete set of 2 k primitive idempotents each with k factors as in (^). 

- The division ring IfC = fC£ Ptq f is isomorphic to R or C or M when (p — q) mod 8 is 
0,1,2, or 3,7 or 4,5,6. 

- The mapping S x K — » S, or (?/>, A) — > ?/>A defines a right Kdinear structure on the 
spinor space S = C£ P:Q f (cf. 0]). 

Example 1. In ~ R(4) we have k = 2 and / = i(l + ei)^(l + e34), = e3e4 = 
e3 A e4 is a primitive idempotent. The ring 1C ~ R is just spanned by {1}r and a real 
basis for S = Cl^xf ma Y be generated by {1, e2, e3, e23}iR (here e23 = e2e3 = e2 A e3.) 
Example 2. In C£ 3fi ~ C(2) we have k = 1 and / = |(1 + ei) is a primitive idempotent. 
The ring K ~ C may be spanned by {l,e23}]R and a basis for S = Cl^flf over IK may be 
generated by {l,e 2 }K- 

Example 3. In Ct\$ ~ H(2), the Clifford polynomial / = |(l+ei4), ei4 = eie4 = ei Ae4, 
is a primitive idempotent. Thus, the ring K ~ EI may be spanned by {1, e2, e3, e23}iR and 
a basis for S = Cl\^f as a right-quaternionic space over K may be generated by {1, 

2 Exponential of a real matrix 

We now proceed to exponentiate a real 4x4 matrix using the spinor representation 7 of 
Ci^i from Example 1. Instead of Clz,x one could also use Ci2,2, the Clifford algebra of 
the neutral signature (2, 2), since (7^2,2 — R(4). From now on = e^e^ = e, A e^, i ^ j, 
K = {Id}w ~ R, and Id denotes the unit element of C£ 3>1 in 'CLIFFORD' 

Recall the following facts about the simple algebra Cl^\ ~ R(4) and its spinor space S: 

4 For the purpose of this paper, it is enough to consider simple Clifford algebras only. 
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C4,i = {l,ej, eij, eijk, e ijk i} K , i < j < k < I, i, j, k,l = 1, . . . ,4. 
5 = C4,i/ = {/i = /, /2 = e 2 /, / 3 = e 3 /, U = e 23 /} K . 

Each basis monomial has a unique matrix 7e i fe! representation in the spinor 
basis fi, i = 1, . . . , 4. For example, the basis 1-vectors ei,e2,e3,e4 are represented 
under 7 as: 





/I 
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-1 
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> 7e 4 = 
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0/ 



(3) 



Since 7 : R(4) — > C?3 1 is a linear isomorphism of algebras, matrices representing Clifford 
monomials of higher ranks are matrix products of matrices shown in (||). For example, 



Tejjfci — 7ei7ej7e fe 7e; : 



7ei234 — 7ei7e2 7e3 7e4 



/ 1 0\ 

-10 

1 

\ -1 0/ 



(4) 



Then, a matrix representing any Clifford polynomial may be found by the linearity of 7. 

Relevant information about Cl^,i is stored in 'CLIFFORD' and can be retrieved as 
follows: 

> restart : with (Cliff 3) :dim:=4:B:=linalg[diag] (1,1,1,-1) : 

> eval (makealiases (dim) ) : data: =clidata() ; 



data :- 



1111 

[real, 4, simple, cmulQ(- Id + - el, — Id + — e34), 

[Id, e2, e3, e23], [Id], [Id, el, e3, e23]] 



In the Maple list data above, 



- real, 4, and simple mean that CI3 \ is a simple algebra isomorphic to M(4). 

- The fourth element data [4] in the list 'data' is a primitive idempotent / written 
as a Clifford product of two Clifford polynomials (Clifford product in orthogonal 
Clifford algebras is realized in 'CLIFFORD' through a procedure 'cmulQ'). 

- The list [Id, el, e3, e23] contains generators of the spinor space S = Cl^^f over the 
reals R (compare with Example 1 above). 
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- The list [Id] contains the only basis element of the field K C CK31, that is, the 
identity element of Cl^x. 

- The final list [Id, e2, e3, e23] contains generators of the spinor space S = C£s t if 
over the field IK. In this case it coincides with data [5] since K ~ R. 

Thus, a real spinor basis in S consists of the following four polynomials: 
> f l:=f ;f2:=cmulQ(e2,f ) ;f3:=cmulQ(e3,f ) ;f4:= cmulQ(e23,f ) ; 

/1:=- Id + - e34 + - el + - el34, /2:=~ e2 + - e234 - - el2 - - el234 
4 4 44 44 4 4 



/3:=- e3 + - e4 - - el3 - - el4, 
4 4 4 4 



/4:=- e23 + - e24 + - el23 + - el24 
4 4 4 4 



(5) 



Procedure 'matKrepr ' allows us now to compute 16 matrices m[i] representing each 
basis monomial in Cl% 1 . 

> for i from 1 to 16 do 

> lprint ('The basis element ', clibas [i] , 

'is represented by the following matrix:); 

> m[i] : =subs (Id=l , matKrepr (clibas [i] ) ) od: 

Let's define a 4 x 4 real matrix A without a complete set of eigenvectors. Therefore, A 
cannot be diagonalized. 

> A :=linalg [matrix] (4,4, [0, 1 ,0,0,-1 ,2,0,0,-1 , 1 , 1 ,0,-1 , 1 ,0, 1] ) ; 

> linalg [eigenvects] (A) ; #A has incomplete set of eigenvectors 



10 
-12 
-1110 
-110 1 



[1, 4, {[0, 0, 1, 0], [1, 1, 0, 0], [0, 0, 0, 1]}] 



(6) 

1 with an algebraic 



Maple output in (j6j) shows that A has only one eigenvalue A 
multiplicity 4 and a geometric multiplicity 3. 

In the Appendix, one can find a procedure 'phi' which gives the isomorphism ip from 
R(4) to Ct^i- It can find the image p = phi(^4) of any real 4x4 matrix A using the 
previously computed matrices m[i]. In particular, the image p of A under <p is computed 
as follows: 



> FBgens : = [Id] ; #assigning a basis element of K 

> p :=phi( A, m, FBgens ) ; #finding the image of A in CI (3,1) 

, 1111111 1 

p:=Id - - el e3 e4 + - el2 e23 e24 - - el34 + - el234 

F 2222222 2 



(7) 



Let's go back to the exponentiation problem. So far we have found a Clifford polynomial 
p in C£s 1 which is the isomorphic image of A. We will now compute a sequence of finite 
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power series expansions of p up to a specified order N. Procedure ' sexp' (defined in 
the Appendix) finds these expansions, which are just Clifford polynomials, modulo the 
minimal polynomial p(x) of p. The minimal polynomial p(x) can be computed using a 
procedure ' climinpoly ' . 

> p(x)=climinpoly(p) ; 

p(x) = x 2 -2x + l (8) 

It can be easily verified that the polynomial (§) is satisfied by p = f(A) and that it is 
also the minimal polynomial of A. 

> cmul(p,p)-2*p+Id; #p satisfies its own minimal polynomial 







> linalg [minpoly] (A,x) ; #matrix A has the same minimal polynomial as p 



x 2 - 2x + 1 



A finite sequence of say 20 Clifford polynomials approximating exp(p) can now be com- 
puted. 

> N:=20:for i from 1 to N do p . i : =sexp(p, i) od:# we want 20 polynomials 

For example, Maple displays polynomial P20 as follows: 

> p_lim:=p.20; 

6613313319248080001 r , 82666416490601 , 82666416490601 

T) iXTfli ' — ± Qj €-1. Co 

1 ' 2432902008176640000 60822550204416 60822550204416 

82666416490601 82666416490601 82666416490601 

~ 60822550204416 6 + 60822550204416 6 ~ 60822550204416 6 

82666416490601 82666416490601 82666416490601 

e24 el34 + el234 

60822550204416 60822550204416 60822550204416 

Having computed the approximation polynomials pi,P2, ■ ■ ■ ,Pn, N = 20, one can show 
that the sequence converges to some limiting polynomictl pum, 

by verifying that \pi~Pj\ < e 
for i,j > M, M sufficiently large, in one of the Maple's built-in polynomial norms. 

Finally, we map back pu m into a 4 x 4 matrix which approximates exp(^4) up to and 
including the terms of order N. 

> expA:=0:for i from 1 to nops(clibas) do 

> expA : =evalm(expA+coef f (p_lim, clibas [i] )*m[i] )od: 

> evalm(expA) ; #the matrix exponent of A 
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2432902008176640000 

-82666416490601 
30411275102208 

-82666416490601 
30411275102208 

-82666416490601 
30411275102208 



82666416490601 
30411275102208 

7775794614048301 
1430277488640000 

82666416490601 
30411275102208 

82666416490601 
30411275102208 






6613313319248080001 









2432902008176640000 

n 6613313319248080001 
u 2432902008176640000 



Although A had an incomplete set of eigenvectors, Maple can find exp(A) in a closed 
form. 

> mA : =linalg [exponential] (A) ; 



mA 



e 

-e 2e 

— e e e 

— e e e 



Notice that our result is very close to the Maple closed-form result: 

> map(evalf , evalm(expA) ) ; 

[.41103176233121648585 10~ 18 , 2.7182818284590452349, 0, 0] 
[-2.7182818284590452349, 5.4365636569180904703, 0, 0] 

[-2.7182818284590452349, 2.7182818284590452349, 2.7182818284590452353, 0] 
[-2.7182818284590452349, 2.7182818284590452349, 0, 2.7182818284590452353] 

The 1-norm of the difference matrix between mA and expA can be computed in Maple as 
follows: 

> evalf (linalg [norm] (mA-expA,l)) ; 

.210~ 17 

3 Exponential of a complex matrix 

In this section we exponentiate a complex 2x2 matrix using a spinor representation of 
C^3,o — C(2) (see Example 2 above). Note that instead of using Cl^fi, one could also use 



Ct\p since Cl\^ ~ C(2). As before, = e^efe = e« A A e&, i,j,k = 1, 



,3, K 



{Id, e 2 3}]R c± C, e| 3 = —Id, where Id denotes the unit element of Ci^fi in 'CLIFFORD'. 
Recall these facts about the simple algebra Ci^jo and its spinor space S: 



CI 



3,0 



{1, e i} e ij} e ijk } R , i < j < k. 



S = C£ 3fi f = {fi = f,f 2 = e 2 f, h = e 3 /, U = e 23 /}i 
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- S = C£ 3fi f = {fi = f,h = e 2 /k- 

For example, the basis 1-vectors are represented in the spinor basis {/i, /2} by these three 
matrices in K(2) well known as the Pauli matrices: 

7ei = (j -l)' 7e2 = (? (">)' 7e3 = (el _ 23 )- (9) 
The following information about Ci^fl is stored in 'CLIFFORD': 

> dim:=3:B:=linalg[diag] (1,1,1) : 

> data: =clidata() ; 

data := [complex, 2, simple, -Id + ~^ e ^i [Id, e2, e3, e23], [Id, e23], [Id, e2]] 

Now we define a Grassmann basis in Ct$ o, assign a primitive idempotent to /, and 
generate a spinor basis for S = C£s t of- 

> clibas : =cbasis (dim) ; #ordered basis in CI (3,0) 

clibas := [Id, el, e2, e3, el2, el3, e23, el23] 

> f:=data[4]; #a primitive idempotent in CI (3,0) 

> sbasis : =minimalideal (clibas ,f , ' lef t ' ) ; #find a real basis in Cl(B)f 

sbasis := 

[[I Id + 1 el, 1 e2 - X - el2, \ ^>~\ d3, ^23 + ^ el23], [Id, e2, e3, e23], Ze/t] 

> f basis : =Kfield(sbasis ,f) ; #find a basis for the field K 

f basis := [[- Id + - el, \ e23 + \ el23 ]i e23 ]] 

> SBgens : =sbasis [2] ; #generators for a real basis in S 

SBgens := [Id, el, e3, e23] 

> FBgens : =f basis [2] ; #generators for K 

FBgens := [Id, e23] 
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In the above, 'sbasis' is a K-basis returned for S = Cl^pf. Since in the current 
signature (3,0) we have K = {Id, e23} R ~ C, cmulQ(e23, e23) = -Id, and C£ 3fi ^ C(2), 
the output from ' spinorKbasis ' shown below has two basis vectors and their generators 
modulo /: 

> Kbasis : =spinorKbasis (SBgens , f ,FBgens , 'left ' ) ; 

Kbasis .- [[I Id + I e l, 1 e 2 - X - e!2], [Id, e2], left] 



> cmulQ(f,f); #verifying that f is an idempotent 

\ Id+ \ el 

Note that the second list in 'Kbasis' contains generators of the first list modulo the 
idempotent /. Thus, the spinor basis in S over IK consists of the following two polynomi- 
als: 

> for i from 1 to nops (Kbasis [1] ) do f . i : =Kbasis [1] [i] od; 



/l:=I/d+iel, /2 :=le2-l e 12 



(10) 



We are in a position now to compute matrices m[i] representing basis elements in C£ 3j q. We 
will only display Clifford-algebra valued matrices representing the 1-vectors {ei,e2,e3} 
and the unit pseudoscalar ei23 = eie2e3. 

> for i from 1 to nops(clibas) do 

> lprint ('The basis element ', clibas [i] , 

'is represented by the following matrix:'); 

> m [i] : =subs (Id=l ,matKrepr (clibas [i] ) ) od: 

The basis element el is represented by the following matrix: 



m 2 :-- 



1 
-1 



The basis element e2 is represented by the following matrix: 



m 3 : = 



1 

1 



The basis element e3 is represented by the following matrix: 
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7714 '■' 



-e23 
e23 



The basis element el23 is represented by the following matrix: 



m 8 :-- 



e23 
e23 



As an example, let's define a complex 2x2 matrix A and let's find its eigenvectors: 

> A :=linalg [matrix] (2,2, [1+2*1,1-3*1,1-1,-2*1] ) ; #defining A 

> linalg [eigenvects] (A) ; 



A :-- 



1 + 21 1-31 
1-1 -21 



[I + iV-23-87, 1, { 



| + IV-23-8/ + /+ij(i + ^V-23-87), 1 



}], 



[\~\ V-23-8/, 1, { 



1 V _ 23 _ 8/ + /+ _ 1 V-23 - 87), 1 



}] 



The image of ^4 in C^o under the isomorphism <£> : C(2) — > C^o can now be computed. 
Recall that ' FBgens ' defined above contained the basis elements of the complex field K 
in C£ 3t0 . 

> evalm(A) ;p:=phi(A,m, FBgens) ; #finding image of A in Cl(3,0) 



1 + 21 1-31 
1-1 -21 



p ■= I Id + i el + e2 + e3 + 2 el3 + 2 e23 

Thus, we have found a Clifford polynomial p in ^3,0 which is the isomorphic image of 
A. We will now compute a sequence of finite power expansions of p up to and including 
power = 30 using the procedure 'sexp'. This sequence of Clifford polynomials should 
converge to a polynomial pu rn , the image under ip of the matrix exponential exp(^4). First, 
we find the real minimal polynomial p(x) of p (called 'pol' in Maple). 

> pol : =climinpoly (p) ; #find the real minimal polynomial of p 



pol ■= x 4 - 2 x 3 + 13 x 1 - 12 x + 40 
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> &c(p$4)-2*&c(p$3)+13*&c(p$2)-12*p+40*Id;#checking that p satisfies pol 



Observe that matrix A has the following complex minimal polynomial 'pol2': 

> pol2 : =linalg [minpoly] (A,x) ; 

pol2 :=6 + 2I -x + x 2 

> evalm(&*(A$2)-A+6+2*I) ; 




Furthermore, since {Id, el23}]R is another copy of the complex field K in Ci^fi, we can 
easily verify that the Clifford polynomial p also satisfies the complex minimal polynomial 
'pol2' of A if we replace 1 with Id and I with el23, namely: 

> &c(p$2)-p+6*Id+2*el23; 



On the other hand, matrix A of course satisfies the polynomial 'pol': 

> evalm(&*(A$4)-2*&*(A$3)+13*&*(A$2)-12*A+40) ; 




As expected, the complex minimal polynomial of A is a factor of the real minimal poly- 
nomial of p: 

> divide (pol, pol2) ; 

true 

> pol3 : =quo (pol , pol2 , x) ; 

pol3:=x 2 -x + 6-21 

Let's check that pol3 * pol2 = pol: 

> pol ; expand (pol3 * pol2) ; 

x 4 - 2x 3 + 13x 2 - 12 x + 40 
x 4 -2x 3 + 13x 2 - 12 x + 40 
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The following loop computes Clifford polynomials pi approximating exp(p) in C£s,o- 
We will only display polynomial p^o and assign it to pu m . 

> Digits: =20: 

> N:=30:for i from 1 to N do p . i : =sexp(p, i) od; 

> p_lim:=p.N: 

739418826545208898275600203389 140606618686769098555631609225939 

P 544108430383981658741145600000 + 176835239874794039090872320000000 6 

13294860446171527820401106221093 ^ 5429376085448859186420447465893 

~ 88417619937397019545436160000000 & + 12631088562485288506490880000000 & 

50830755859220399836279191881837 15796535483801410769637551225479 

+ 44208809968698509772718080000000 & + 22104404984349254886359040000000 6 

537129223345642211370021843709 24569201649575451209456052913 

■ el23 el2 



1184164552732995797483520000000 84691206836587183472640000000 

By picking up numeric coefficients of the basis monomials in the subsequent approxima- 
tions to exp(p), one can get an idea about the approximation errors. 

> sort( [op(L:=cliterms(p_lim))] ,bygrade) : 

> for i from 1 to nops(L) do 

> L. i :=map(evalf , [seq(coef f (p. j ,L [i] ) , j= 1..N)]) od: 

> approxerror : = 

max(seq(min(seq(abs(L. j [i]-L. j [i-1] ) , i=2..N)), j=l . .nops (L) ) ) ; 



approxerror := .1 10 



-19 



Having computed the finite sequence of polynomials pi one can again show by using 
Maple's built-in polynomial norm functions that this is a convergent sequence. For exam- 
ple, in the infinity norm one gets |p29 — P3o\ < -6 x 10~ 20 and \pi — Pj\ — > as i,j — ► oo. 

Thus, we have found an approximation pn m to the power series expansion of exp(p) 
in Cisfi up to and including terms of degree iV = 30. Finally, we map back pn m into 
a 2 x 2 complex matrix which approximates exp(A). We expand pu m over the matrices 
m[i]: 

> expA:=0: for i from 1 to nops(clibas) do 

> expA:=evalm(expA+coeff (p_lim, clibas [i] )*m[i]) od: 

> evalm(expA) ; #the matrix exponent of A 

7121749995744556670281318348249 596909308415457533523842428577 



12631088562485288506490880000000 2286662584587853953761280000000 
7789021393665659776614645092453 6228189267183180110479443301 



17683523987479403909087232000000 394281471292740332421 1200000 
12355386075985243242271013020079 340405770696784948489921131029 



88417619937397019545436160000000 472821496991427912007680000000 
31743144776163499207933472943947 5959141766058476626587221301857 



14736269989566169924239360000000 5101016534849828050698240000000 



e23 
e23, 
e23 
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Maple can find the exponent of A in a closed form with its ' linalg [exponential] ' 
command. We won't display the result but we will just compare it numerically with our 
result saved in ' expA ' . 

> mA :=linalg [exponential] (A) : 

Let's replace the monomial e2we3 in 'expA' with the imaginary unit 7 used by Maple 
and let's apply 'evalf ' to the entries of ' expA': 

> f expA : =subs (e2we3=I , map (evalf , evalm(expA) ) ) ; 

fexpA : = 

[-.56382709696901085353 + .26103952215715461164 / , 

- .44046771442052942162 - 1.5796302186766888057/] 
[.13973895796712599250 - .71994563035478140661/, 

- 2.1540827359052753813 - 1.1682263182928324795 7] 

> fmA : =map( evalf ,mA) ; #applying 'evalf' to mA 

fmA : = 

[-.56382709696901085362 + .26103952215715461158 1 , 

- .44046771442052942180 - 1.5796302186766888058 7] 
[.13973895796712599243 - .71994563035478140663 7 , 

- 2.1540827359052753816 - 1.1682263182928324795 7] 

Let's check the 1-norm of the difference matrix between 'fmA' and 'fexpA': 

> evalf (linalg [norm] (f mA-f expA , 1) ) ; 

.5059126028 10~ 18 

The floating-point approximation 'fexpA' to exp(yl) is within approximately .5 x 10 -18 
in the matrix || • ||i norm to the closed matrix exponential computed by Maple. 

4 Exponential of a quaternionic matrix 

In order to exponentiate a quaternionic 2x2 matrix, we will use the spinor representation 
of C£i 5 3 ~ H(2) (see Example 3 above). Note that two other algebras could be used 
instead of C£\ :3 , namely, C£o,4 and Cl±$ since both are isomorphic to H(2). As before 
eij = e^ej = A ej, i,j = 1, . . . , 4, but this time K = {Id, e3, 623}^ — H. 

Recall the following facts about the simple algebra Ci\^ and its spinor space S: 

- C£ lj3 = {l,ei,eij,eijk,eijki}igL, i < j < k < I. 

- S = C£ h3 f = {h = f,f2 = e 2 f, h = e 3 /, h = e 2 3/"k. 



- S = C£ li3 f = {f 1 = f,f 2 = e 1 f} K . 
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For example, the basis 1-vectors ei,e2,e3,e4 are represented by: 

7ei = (? o)' 7e2 = (o -e 2 )' 7e « = (o 3 -e 3 )' 7e4= (? V)- (U) 
In order to compute the spinor representation of Ci\^, we proceed as follows: 

> data:=clidata(linalg[diag] (1,-1,-1,-1)) ; 

1 1 

data := [quaternionic, 2, simple, — Id + ~^ e ^i [Id, el, e2, e3, el2, el3, e23, el23], 
[Id, e2, e3, e23], [Id, el]] 

We define a Grassmann basis in Ct\$, assign a primitive idempotent to /, and generate a 
spinor basis for S = Cl\$j ' . 

> clibas : =cbasis (dim) ; #ordered basis in Cl(l,3) 

clibas := 

[Id, el, e2, e3, e4, el2, el3, el4, e23, e24, e34, el23, el24, el34, e234, el234] 

> f:=data[4]; #a primitive idempotent in Cl(l,3) 

/:=5«+iel4 

Next, we compute a real basis in the spinor space S = Cii^f using the command 
' minimalideal ' : 

> sbasis : =minimalideal (clibas ,f , ' lef t ' ) ;#f ind a real basis in Cl(B)f 

sbasis := \\- Id+- el4, - el + - e4, - e2 - - el24, - e3 - - el34, - el2 - - e24, 
LL 2 2 ' 2 2'2 2 2 2 2 2' 

- el3 - - e34, - e23 + - el234, - el23 + - e234], 
2 2 ' 2 2 '2 2 J 

[Id, el, e2, e3, el2, el3, e23, el23], left] 

In the following, we compute a basis for the subalgebra K: 

> f basis : =Kfield(sbasis ,f) ; #a basis for the field K 

fbasis := 

[[I Id+\ el4, \e2-\ el24, \ e3 - \ el34, \ e23 + \ el234], [Id, el, e3, e23]] 



> SBgens : =sbasis [2] ; #generators for a real basis in S 

SBgens := [Id, el, e2, e3, el2, el3, e23, el23] 
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Thus, a possible set of generators for K is: 
> FBgens : =f basis [2] ; #generators for K 

FBgens := [Id, e2, e3, e23] 



(12) 



In the above, ' sbasis' is a real basis for S = C£\^f. Since in the current signature 
(1,3) we have that IK = {Id, e2, e3, e23} R ~ M and C£ lj3 = 11(2), the output from 
' spinorKbasis ' shown below has two basis vectors and their generators modulo / for S 
over IK: 

> Kbasis : =spinorKbasis (SBgens ,f , FBgens , ' left ' ) ; 

Kbasis := [[- Id+- el4, - el + - e4], [Id, el], Ze/t] 

> cmulQ(f,f); #f is an idempotent in Cl(l,3) 

\ld + \elA 

Notice that the generators of the first list in 'Kbasis' are listed in Kbasis [2]. Further- 
more, a spinor basis in S over IK consists of the following two polynomials f\ and f 2 '- 

> for i from 1 to nops (Kbasis [1] ) do f . i : =Kbasis [1] [i] od; 



n:=\ld+\e\A, / 2:=i e l + i e 4 



(13) 



Using the procedure 'matKrepr' we can now find matrices m[i] with entries in IK rep- 
resenting basis monomials in Cl\^. Below we will display only matrices representing the 
1-vectors ei,e2,e3 and e^. 

> for i from 1 to nops(clibas) do 

> lprint ('The basis element clibas [i] , 

'is represented by the following matrix:'); 

> m [i] : =subs (Id=l , matKrepr (clibas [i] ) ) od; 

The basis element el is represented by the following matrix: 



m 2 :- 



1 

1 



The basis element e2 is represented by the following matrix: 



m 3 



e2 
-e2 
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The basis element e3 is represented by the following matrix: 



7B4 : = 



e3 
-e3 



The basis element e4 is represented by the following matrix: 



m 5 :- 



-1 

1 



Let's define a 2 x 2 quaternionic matrix A. In Maple, we will represent the standard 
quaternionic basis {l,i,j,k} as {1, 'ii' , ' jj ' , 'kk'}. Later we will make substitutions: 
' ii' — * e2, 'jj' — ► e3, 'kk' — > e2we3 since, as we may recall from Example 3 above, 
IK = {l,e2,e 3 ,e 2 3} R . 

> A :=linalg [matrix] (2,2, [l+2*'ii'-3*'kk' ,2+'ii' -2*'jj', 

> 'kk'-3*'ii' ,2*'kk'-2*' j j '] ) ; #defining a quaternionic matrix A 



A := 



1 + 2 ii - 3 kk 2 + ii-2jj 
kk — 3ii 2kk — 2jj 



(14) 



The isomorphism if : H(2) — ► Cl\$ has been defined in Maple through the procedure 
'phi ' (see the Appendix). This way we can find image p in Ci\ 3 of any matrix A. Recall 



that 'FBgens' in (12) contains the basis elements of the field IK. 

> p:=phi(A, m, FBgens) ;#finding image of A in Cl(l,3) 

p ■= - Id + el + e2 + e3 - e4 - 2 el2 + el3 + - el4 - ^ e23 + e24 + e34+ 

1 15 

- el23 - el24 + el34 + - e234 - - el234 

2 2 2 

The minimal polynomial p(x) of p in Ct\$ is then found with the procedure ' climinpoly ' : 

> climinpoly (p) ; 

x 4 - 2x 3 + 16x 2 + 10 x + 330 



So far we have found a Clifford polynomial p in Cl\$ which is the isomorphic image 
of the quaternionic matrix A. We will now compute a sequence of finite power expansions 
of p using the procedure 'sexp'. This sequence of Clifford polynomials will be shown to 
converge to a polynomial pum that is the image of exp(A). For example, polynomial p20 
= sexp (p, 20) looks as follows: 

> for i from 1 to 20 do p . i : =sexp(p , i) od; 
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+ 
+ 
+ 



68240889697169513 50515123107772493 

P '~ ~ 10861169679360000 9503523469440000 6 

976049744897473 76665127748453 

el23 e234 

638892334080000 66691392768000 

23336382714907219 1736342897976643 

el24 el34 

152056375511040000 1974758123520000 

9030311044661089 t 802551523836832291 

el234 H el2 

1407929402880000 152056375511040000 

907882088300711 4304638284278411 

el3 H e23 

365520133440000 4472246338560000 

360072975386539 ^ 19812017405738017 

~ 116162242560000 & ~ 76028187755520000 6 

1889118161676113 277471312336316837 

~ 703964701440000 &l ~ 152056375511040000 62 

98120514192871531 25277099300039 

~ 152056375511040000 ^ ~ 44722463385600 6 



Thus, we have a finite sequence of Clifford polynomials pi approximating exp(p). Next, 
for each of the 16 basis monomials present in all polynomials, we create a sequence Sj (or 
sj in Maple) of its coefficients. 

> for j from 1 to nops(clibas) do 

> s . j : =map(evalf , [seq(coef f (p . i , clibas [j] ) , i=l . . N)] ) od: 

For example, the sequence si of the coefficients of the identity element Id is: 

> si; 



[1.500000000, -2., -6.916666667, -18.66666667, -20.22500000, -10.85972222, 
-5.099206349, -3.980456349, -5.027722663, -6.129274691, -6.428549232, 
-6.368049418, -6.301487892, -6.280796253, -6.280315663, -6.282290205, 
-6.282986035, -6.283054064, -6.283026981, -6.283014787] 



Having computed the finite sequence of polynomials pi,P2, ■ ■ ■ ,P20, one can again verify 
that this is a convergent sequence by using any of the Maple's built-in polynomial norm 
functions to estimate norms of the differences pi — pj for i,j = 1, . . . , 20. It can be again 
observed that \pi —pj | — > as i, j — > oo. Finally, we map back pum — P20 into a 2 x 2 matrix 
'expA' which approximates exp(^4) up to and including terms of order N = 20. After 
expressing back the basis elements {Id, e2, e3, e2we3} in terms of {1, ' ii' , ' jj ' , 'kk'} 
we obtain: 

> p_lim:=p20: 

> expA:=0:for i from 1 to nops (clibas) do 

> expA : =evalm(expA+coef f (p_limit , clibas [i] ) *m [i] ) od: 

> sexpA:=subs({e2we3='kk' ,e3=' jj ' ,e2='ii'}, evalm(expA) ) ; 



sexpA := 
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58889470322671 301630543173.. 1778894447566499.. 
" 8999548740000 ~ 152472320000 %% + 7602818775552000 33 
560815647244431793 , , 10065855790684619 5520266650930879 . . 
76028187755520000 ' 4751761734720000 2534272925184000 
748687448521121 .. 203548165276035707 , ,' 
95995186560000 JJ 76028187755520000 

30874478783885813 33523343384679259.. 3844312687422001.. 
" 9503523469440000 + 4001483566080000 " + 1357646209920000 33 
2613788546323897 , , 228937105237224287 127067464810704809 
6911653432320000 ' 38014093877760000 76028187755520000 
38636486222845507 . . 414457945578965819 
+ 25342729251840000 33 ~ 76028187755520000 

> f expA : =map(evalf , evalm(sexpA) ) ; #f loating-point approximation 

fexpA : = 

[-6.543602577 - 1.978264272 ii + .2339782783 jj + 7.376417403 kk , 
-2.118341860 - 2.178244733 + 7.799218642 jj + 2.677272355 kk] 
[-3.248740205 + 8.377728618 ii + 2.831601237 jj + .3781712396 kk , 
-6.022426997 + 1.671320448 ii + 1.524559010 jj - 5.451372153 kk] 

Thus, matrix ' sexpA ' is the exponential of the quaternionic matrix A from (|l^ ) computed 
with the Clifford algebra Cl\ 3. 



5 Conclusions 

We have translated the problem of matrix exponentiation e A , A 6 K(n), into the problem 
of computing e p in the Clifford algebra C£(Q) isomorphic to K(n). This approach, alter- 
native to the standard linear algebra methods, is based on the spinor representation of 
C£(Q). It should be equally applicable to other functions representable as power series. 
Another use for the isomorphism between C£ Pjq and appropriate matrix rings could be to 
finding the Jordan canonical form of A in terms of idempotent and nilpotent Clifford poly- 
nomials from C£(Q) (see also flfj]] and 0] for 

more on the Jordan form and its relation to 
the Clifford algebra). Generally speaking, any linear algebra property of A can be related 
to a corresponding property of p, its isomorphic image in C£(Q), and it can be stated in 
the purely symbolic non-matrix language of the Clifford algebra. These investigations are 
greatly facilitated with 'CLIFFORD'. At pj interested Reader my find complete Maple 
worksheets with the above and other computations. 
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7 Appendix 

The procedures described in this Appendix will work provided the Maple package 'CLIF- 
FORD' has been loaded first into a worksheet.^] Procedure 'phi' was used above to 
provide the isomorphism tp between the matrix algebras M(4), C(2), and H(2) and, respec- 
tively, the Clifford algebras C-^i, C^o, and Ci\$. 

> phi :=proc(A: :matrix,m: : table, FBgens : : list (climon) ) 
local N, n, cb, f b, AA, M, a, j ,L,sys,vars , sol, p; global B; 
if nops (FBgens) =1 then AA:=evalm(A) elif 

nops (FBgens) =2 then f b : =op (remove (has , FBgens , Id) ) ; 
AA:=subs(I=fb,evalm(A)) elif 

nops (FBgens) =4 then fb:=sort (remove (has, FBgens, Id) ,bygrade) ; 

AA:=subs('ii'=fb[l] , ' jj '=fb[2] , 'kk'=fb[3] ,evalm(A)) 

else ERROR('wrong number of elements 'FBgens'') fi; 
N:=nops( [indices (m)] ) ;n:=linalg[coldim] (B) : cb:=cbasis(n) ; 
M:=map(displayid,evalm(AA-add(a[j] *m[j] , j=l. .N))) ; 
L : =map(clicollect , convert (M,mlist) ) ; 

sys : =op(map(coef f s ,L, FBgens) ) ; vars :=seq(a[j] ,j=l. . N) ; 
sol : =solve ( sys , vars) ; vars : =seq(a [j] *cb [j] , j=l . . N) ; 
p:=subs(sol,p) ; RETURN (p) 
end: 

Procedure 'climinpoly' finds a real minimal polynomial of any Clifford polynomial p in 
an arbitrary Clifford algebra Ci P:q . 

> climinpoly :=proc(p: : clipolynom, s : : string) 

local dp , L , f lag , pp , expr , a , k , eq, sys , vars , sol , poly ; 
option remember; 

dp : =displayid(p) :L:=[Id,dp] ;flag:=false: 
while not flag do 

pp : =cmul (L [nops (L) ] , dp) : 

expr : =expand(add(a[k] *L [k] ,k=l . .nops (L) ) ) ; 

eq: =clicollect (pp-expr) ; sys : =coef f s (eq, cliterms (eq) ) ; 

vars : =seq(a[k] ,k=l . .nops (L) ) ; sol : =solve (sys , vars) : 

if solo then flag:=true else L: = [op(L) ,pp] fi; 

od; 

poly :='x' "nops (L)-add(a[k] *'x' " (k-1) ,k=l . .nops(L) ) ; 
if nargs=l then RETURN (sort (subs (op (sol) , poly))) 

else RETURN( [sort (subs (op(sol) ,poly)) ,L] ) fi; 
end: 

n 

Procedure 'sexp' finds a finite formal power series expansion ^ (p k /k\) of any Clifford 

k=0 

polynomial p up to and including the degree specified as its second argument. Com- 
putation of the powers of p in C£ Ptq is performed modulo the real minimal polynomial 
of p. 

> sexp : =proc (p : : clipolynom, n: :posint) local i,d,L,Lp,pol,poly,k; 
pol : =climinpoly (p , 's') ;readlib(powmod) ; 
poly:=add(powmod( 'x' ,k,pol [1] , 'x' )/k! ,k=0. .n) ; 

L : = [op (poly) ] ; Lp : = [] : 

for i from 1 to nops(L) do 



5 To download 'CLIFFORD', see the Web site in §. 
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d:=degree(L[i] ) ; 

if d=0 then Lp : = [op(Lp) ,L [i] *Id] else 

Lp:=[op(Lp) ,coeffs(L[i])*pol[2] [d+1]] fi od; 
RETURN ( add ( Lp [ i ] ,i=l. .nops(Lp))) 
end: 
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